---
  title: "CARRS RDD analysis :Robustness check keeping only follow uo observations within 12 months"
knit: (function(input, ...) {rmarkdown::render(input, output_file="")})
output: pdf_document
geometry: "left=2cm,right=2cm,top=1cm,bottom=1cm"
date: "`r format(Sys.time(), '%d %B, %Y')`"
---
  
```{r , setup, include=FALSE, eval = TRUE}

knitr::opts_chunk$set(
  comment = '', warning=FALSE,  message=FALSE, results= FALSE , fig.width = 9,fig.height = 7, echo=FALSE, root.dir =""
)

knitr::opts_knit$set(
 root.dir = ""    
)
```


```{r packages}
rm(list=ls())
library(tidyverse)
library(gridExtra)
library(dplyr)
library(rdrobust)
library(rddensity)
library(viridisLite)
library(foreign)
library(ggplot2)
library(data.table)
library(doBy)
library(socviz)
library(openxlsx) 
library(modelsummary)
library(labelled)
library(gtsummary)
library(gt)
library(kableExtra)
```

```{r loaddata}
data_treatment<-data.frame(read.csv("carrs_full.csv"))
data_treatment<-as.data.table(data_treatment)

data_diagnosis<-data.frame(read.csv("carrs_selected.csv"))
data_diagnosis<-as.data.table(data_diagnosis)

#sdiff calculates the number of months between the two surveys
data_diagnosis<-data_diagnosis[sdiff>=0 & sdiff<13, ]
data_treatment<-data_treatment[sdiff>=0 & sdiff<13, ]
```
##calculating the bandwidth for each sub group

```{r}
# Males and females separately
rv.update = function(b,c) {
    out=c[female==b,runningvar:=ifelse((((sys-140)/(sd(sys))) == ((dias-90)/(sd(dias)))),
                                            (sys-140)/(sd(sys)),
                                            ifelse((((sys-140)/(sd(sys))) > ((dias-90)/(sd(dias)))),((sys-140)/(sd(sys))),((dias-90)/(sd(dias)))))]
  a_b<<-out
}

rv.update(1,data_diagnosis)
diag_female<-a_b[female==1,]
rm(a_b)

rv.update(0,data_diagnosis)
diag_male<-a_b[female==0,]
rm(a_b)

rv.update(1,data_treatment)
treat_female<-a_b[female==1,]
rm(a_b)

rv.update(0,data_diagnosis)
treat_male<-a_b[female==0,]
rm(a_b)

# total population
rv.updateMT = function(c) {
  out=c[,runningvar:=ifelse((((sys-140)/(sd(sys))) == ((dias-90)/(sd(dias)))),
                                             (sys-140)/(sd(sys)),
                                            ifelse((((sys-140)/(sd(sys))) > ((dias-90)/(sd(dias)))),((sys-140)/(sd(sys))),((dias-90)/(sd(dias)))))]
  a_b<<-out
}

rv.updateMT(data_diagnosis)
diag_total<-a_b
rm(a_b)

rv.updateMT(data_treatment)
treat_total<-a_b
rm(a_b)
```
##RDD Outputs function

```{r fuc}
carrs_rdd = function(a,b,d) {
  d$Y<-a
  d$X<-b
  
  out =  rdrobust(y = d$Y, d$X, kernel = "triangular", p = 1, covs = cbind(d$w01, d$w23, d$w45), 
                  c = 0, cluster = d$pid, rho = 1 , all = T)
  summary(out)
  
  ret <- data.frame(
    sex = c(as.character(round(100*out$coef[1,1],2)), as.character(round(out$pv[3,1],2)), 
            paste0( "(", as.character(round(100*out$ci[3,1],2)), " - ", as.character(round(100*out$ci[3,2],2)), ")" ), 
            as.character(round(out$bws[1,1],2)), as.character(out$N_h[1] + out$N_h[2]))
  )
  rownames(ret) <- c("Coefficient", "p-value","95\\% CI", "Bandwidth", "Obs. within bandwidth")
  
  output<<-ret %>% add_rownames  
}
```

##Calling the specific subroups and saving the output


```{r diagnosis}
rdd = function(a) {
  carrs_rdd(a$hypdiag,a$runningvar,a)
}

rdd(diag_male)
diag_male <-output%>%mutate(Male=sex)%>%select(-sex)

rdd(diag_female)
diag_female <-output%>%mutate(Female=sex)%>%select(-sex)

rdd(diag_total)
diag_total <-output%>%mutate(Total=sex)%>%select(-sex)

```

```{r treatment}
rdd = function(a) {
  carrs_rdd(a$hypdrug,a$runningvar,a)
}

rdd(treat_male)
treat_male <-output%>%mutate(Male=sex)%>%select(-sex)

rdd(treat_female)
treat_female <-output%>%mutate(Female=sex)%>%select(-sex)

rdd(treat_total)
treat_total <-output%>%mutate(Total=sex)%>%select(-sex)

```


```{r}
total_table<-merge(diag_total,treat_total, by="rowname")
total_table<-total_table %>% arrange(match(rowname, c("Coefficient", "p-value","95\\% CI", "Bandwidth", "Obs. within bandwidth")))

colnames(total_table)[2] <- "Diagnosis"   
colnames(total_table)[3] <- "Treatment"   

female_table<-merge(diag_female,treat_female, by="rowname")
female_table<-female_table%>% arrange(match(rowname, c("Coefficient", "p-value","95\\% CI", "Bandwidth", "Obs. within bandwidth")))

colnames(female_table)[2] <- "Diagnosis"   
colnames(female_table)[3] <- "Treatment"    

male_table<-merge(diag_male,treat_male, by="rowname")
male_table<-male_table %>% arrange(match(rowname, c("Coefficient", "p-value","95\\% CI", "Bandwidth", "Obs. within bandwidth")))

colnames(male_table)[2] <- "Diagnosis" 
colnames(male_table)[3] <- "Treatment"    
```

##table making function

```{r}
  filename <- paste0("limitsurvey_combined.tex")
  caption=paste0("Effect of home-based screening restricted to follow-up surveys conducted within the 12 month recall period")

comb_table<-rbind(total_table,female_table,male_table)

combined<-kbl(comb_table,caption=caption, booktabs = T, linesep = '', escape = FALSE,align = "lcccc", format = "latex",col.names = c(" ","Diagnosis","Treatment"),digits = 2) %>%
    kable_styling(latex_options = "HOLD_position", position = "center") %>%
    pack_rows("Panel A: Total", 1,5) %>%
    pack_rows("Estimation results", 1,3) %>%
    pack_rows("Bandwidth details", 4, 5) %>%
    pack_rows(" ", 6,10) %>%
    pack_rows("Panel B: Female", 6,10) %>%
    pack_rows("Estimation results", 6, 8) %>%
    pack_rows("Bandwidth details", 9, 10) %>%
    pack_rows(" ", 11,15) %>%
    pack_rows("Panel C: Male", 11,15) %>%
    pack_rows("Estimation results", 11, 13) %>%
    pack_rows("Bandwidth details", 14, 15) %>%
    footnote(general = c("The model estimation included a local linear trend, triangular Kernel weights, and a data-driven bandwidth selection (mean squared error optimal bandwidth). The table displays the local average treatment effect in percentage points for individuals with a standardized blood pressure of zero and CIs resulting from robust bias-corrected standard errors clustered at the individual-level.","Abbreviation: CI = confidence interval"),general_title="",threeparttable = T, fixed_small_size = T)%>%
    save_kable(filename, format = "latex")
```


